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This  paper  conducts  a  three-dimensional  (3D)  modeling  study  to  investigate  the  hydrogen  absorption 
process  and  associated  mass  and  heat  transport  in  a  metal  hydride  (LaNi5)  hydrogen  storage  tank.  The 
3D  model  is  further  implemented  numerically  for  validation  purpose  and  the  detailed  investigation  on 
absorption  process.  Results  indicate  that  at  the  very  initial  absorption  stage  the  bed  temperature  evolves 
almost  uniformly,  while  it  varies  greatly  spatially  at  the  latter  stage.  At  the  initial  seconds,  most  hydrogen 
is  absorbed  in  the  region  near  the  cooling  wall  due  to  the  better  heat  removal.  The  absorption  in  the  core  is 
slow  at  the  beginning,  but  becomes  important  at  the  very  end  stage.  It  also  shows  that  the  initial  hydrogen 
flow  in  the  bed  is  several-fold  larger  than  the  latter  stage  and  the  flow  may  provide  extra  cooling  to  the 
hydriding  process.  By  analyzing  the  Peclet  number,  we  find  that  the  heat  convection  by  the  hydrogen  flow 
may  play  an  important  role  in  local  heat  transfer.  This  work  provides  an  important  platform  beneficial  to 
the  fundamental  understanding  of  multi-physics  coupling  phenomena  during  hydrogen  absorption  and 
the  development  of  on-board  hydrogen  storage  technology. 

©  2009  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Hydrogen  as  an  alternative  energy  source  has  been  receiving 
great  attention  worldwide  due  to  the  growing  demand  for  environ¬ 
mentally  clean  fuels  to  replace  nonrenewable  carbonaceous  fuels 
such  as  gasoline  [  1  ].  Hydrogen  storage  is  one  of  the  central  technical 
barriers  to  the  wide  deployment  of  hydrogen  energy,  particularly 
for  automotive  applications  (e.g.  hydrogen  fuel  cell  vehicles).  Cur¬ 
rent  hydrogen  storage  methods  including  gas  compression  and 
liquefaction  are  not  optimal  for  onboard  applications  due  to  the 
associated  large  quantity  of  energy  consumption  and  safety  issues. 
A  promising  alternative  is  solid-metal  hydride  storage:  utilizing 
metal  hydrides  to  absorb/desorb  hydrogen  at  relatively  low  pres¬ 
sures  (<20bar),  which  offers  safety  and  cost  advantage,  however, 
the  obvious  disadvantage  of  weight  for  hydrogen  onboard  stor¬ 
age. 

Developing  effective  solid-state  hydrogen  storage  for  trans¬ 
portation  is  crucial  for  the  success  of  the  hydrogen  economy. 
Besides  energy  density  and  system  cost  requirements,  one  impor¬ 
tant  requirement  by  transportation  applications  is  fast  kinetics  of 
hydrogen  storage,  i.e.  quick  charge  and  discharge.  For  instance,  refu¬ 
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eling  time  for  automobiles  must  be  short  in  minutes  in  order  to 
compete  with  traditional  combustors.  In  addition,  dynamic  charac¬ 
teristics  of  hydrogen  release  should  satisfy  the  need  of  the  fast  load 
variation.  To  develop  such  effective  storage  systems,  physics  of  the 
transport  process  coupled  with  reaction  kinetics  are  highly  in  need, 
such  as  hydrogen  mass  flow  in  the  hydride  bed,  heat  transfer  within 
the  bed,  and  local  hydrogen  absorption  rate. 

Fundamental  modeling  of  hydrogen  storage  has  been  attempted 
by  several  groups.  Jemni  and  Nasrallah  [2]  developed  a  two- 
dimensional  model  for  the  hydrogen  storage  LaNi5-H2  system. 
Efforts  were  also  made  to  experimentally  determine  the  reaction 
kinetics,  equilibrium  pressure,  and  thermal  conductivity  [3].  Naka- 
gawa  et  al.  [4]  present  a  two-dimensional  model  for  heat  and  mass 
transfer  within  the  LaNi5  alloy  bed.  They  found  that  the  absorp¬ 
tion  rate  may  be  limited  by  heat  transfer  due  to  the  low  thermal 
conductivity  of  the  alloy  bed,  therefore  external  cooling/heating 
may  play  an  important  role  in  the  hydriding/dehydriding  process. 
Mat  and  Kaplan  [5]  developed  a  two-dimensional  model  consist¬ 
ing  of  heat  and  mass  transfer  and  chemical  reaction  in  the  LaNi5 
storage  tank  bed.  They  found  that  the  hydride  formation  enhances 
at  lower  equilibrium  pressure  regions  and  the  hydrogen  absorp¬ 
tion  rate  slows  down  due  to  the  reaction  heat  that  increases  the 
bed  temperature.  Other  modeling  approaches  include  Refs.  [6-9], 
which  investigate  the  hydriding/dehydriding  process  and  examine 
the  factors  that  affect  the  tank  performance. 
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Nomenclature 

Bi  the  Biot  number 

cp  specific  heat,  J  kg-1 I<-1 

D  diameter  of  molecule  (m) 

Ea  activation  energy 

H/M  hydrogen  to  metal  atomic  ratio 

AH°  reaction  heat  of  formation  (J  kg-1 ) 
h  heat  transfer  coefficient  (W  m-2 I<-1 ) 

k  thermal  conductivity  (W  nrr 1  K-1 ) 

M  molecular  weight  (kg  mol-1 ) 

n  the  direction  normal  to  the  surface 

P  pressure  (bar) 

Pr  Prandtl  number 

r  radius  (m) 

R  universal  gas  constant,  8.134J  mol-1 K-1 

Re  Reynolds  number 

5  source  term 

t  time  (s) 

T  temperature  (K) 

u  velocity  vector  (ms-1) 

Greek  letters 

a  thermal  diffusivity  (m2  s-1 ) 

p  density  (kg  m-3) 

s  porosity 

Superscripts  and  subscripts 
eq  equilibrium 

eff  effective  value 

g  gas  phase 

m  mass 

s  solid;  saturation 


Though  great  efforts  have  been  made  in  the  hydrogen  storage  of 
metal  hydride  beds,  study  of  detailed  distribution  and  evolution  is 
highly  in  need  to  investigate  the  local  hydriding  process,  coupled 
transport  and  associated  limiting  factors  as  well  as  the  characteris¬ 
tics  of  different  hydriding  stages.  This  paper  seeks  to  develop  a  3D 
model  of  the  hydriding  process  coupled  with  mass/heat  transport 
in  the  hydrogen  tank  and  explore  these  issues  that  are  extremely 
important  for  efficient  hydrogen  storage.  The  model  is  validated 
with  published  experimental  data  and  further  applied  to  investi¬ 
gate  mass/heat  transport  and  reaction  kinetics  during  the  hydriding 
process. 


2.1.  Mass  transport 


Gaseous  hydrogen  flows  through  the  bed  and  reacts  with  the 
hydride  metal.  The  mass  conservation  of  hydrogen  can  be  written 
as 

^  +  V  •  ftu  =  -Sm  (1) 

where  Sm  is  the  local  hydrogen  absorption  reaction  rate  and  u  the 
superficial  velocity.  The  metal  hydride  tank  is  typically  character¬ 
ized  as  a  porous  medium  packed  with  alloy  particles.  The  Darcy’s 
law  can  be  applied  to  describe  the  hydrogen  flow: 

Pgu  =  -§(VP-pgg)  (2) 

where  the  permeability  I<  is  calculated  through  the  Blake-Kozeny 
equation  that  relates  I<  to  the  alloy  particle  radius  and  bed  porosity 
[10]: 


Dj  s 


150(i_£): 


(3) 


where  Dp  is  the  particle  diameter  and  e  the  porosity.  For  a  typical 
size  of  particles  (Dp  ~  1  mm)  and  £  ~  0.5,  I<  is  10-8  to  10-9  m2  s-1 . 

In  the  expansion  region  where  there  are  no  solid  particles,  the 
Navier-Strokes  equation  applies.  A  single-domain  description  can 
be  formulated  with  a  general  moment  equation: 


1 
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dpgU 

~df~ 


V  • 


-VP  +  V  •  r  +  SU  +  psg 


(4) 


where  u  is  the  superficial  velocity.  The  source  term  Su  describes  the 
Darcy’s  forces  (-(/jl/K)u),  which  is  zero  in  the  expansion  region. 
Note  that  Eq.  (1)  gives  a  general  expression  for  the  mass  conserva¬ 
tion  that  can  be  applied  in  both  porous  bed  and  expansion  region, 
e.g.  for  the  expansion  region,  the  porosity  is  1  and  the  source  term 
5m  becomes  zero.  The  gas  density  is  calculated  through  the  ideal 
gas  law: 


pS  _ 

F  RT 


(5) 


2.2.  Heat  transport 


The  energy  equation  in  the  porous  bed  can  be  expressed  as 

^l+v.  (pSc$uT)  =v.(k^vr)  +  sT  (6) 

where 

pcp  =  (1  -  £)pscsp  +  £pscf),  keff  =  (1  -  £)tsks  +  £zsJ<g,  and 
ST  =Sm  [AH°-(cf-cp]  (7) 


2.  Mathematical  modeling 

During  the  hydriding  process,  hydrogen  is  injected  to  a  stor¬ 
age  tank  and  transported  via  the  interstitial  space  of  the  packed 
alloy  particles,  reacting  with  the  unconsolidated  metal  particles. 
In  addition  to  mass  transport,  the  hydriding  process  is  exother¬ 
mic  and  temperature  has  profound  impacts  on  reaction  kinetics, 
therefore  heat  production  and  transfer  within  the  tank  must  be 
included  to  investigate  the  hydrogen  adsorption.  To  develop  a  math¬ 
ematical  model  that  describes  the  hydriding  process,  the  following 
assumptions  are  made:  (1 )  the  flow  bed  can  be  treated  as  a  homoge¬ 
neous  porous  media  with  uniform  morphological  properties,  such 
as  porosity,  tortuosity  and  permeability;  (2)  the  ideal  gas  law  holds 
true  in  the  gas  phase;  (3)  local  solid  and  fluid  are  in  thermal  equi¬ 
librium. 


Note  that  the  thermal  conductivity  of  the  solid  is  much  higher 
than  hydrogen,  therefore  the  contribution  from  hydrogen  can  be 
neglected.  In  addition,  the  interfacial  contact  resistance  among  the 
metal  particles  may  contribute  to  a  considerable  portion  of  the 
thermal  resistance,  which  is  included  in  the  tortuosity  of  the  solid 
matrix  rs.  In  the  literature,  the  value  of  keff  is  around  1 W m-1  K-1. 

2.3.  Hydriding  reaction  kinetics 

The  following  reaction  takes  place  at  the  absorption  sites  of  the 
alloy  (LaNi5): 

LaNi5  +  3.25H2  LaNi5H6  5  (8) 

As  the  alloy  particles  are  solid,  the  gaseous  hydrogen  first  enters 
the  particle  via  the  surface.  The  H  atom  within  the  metal  alloy 
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may  vibrate  around  alloy  atom  forming  H-M  bond  or  take  multiple 
hops  among  metal  alloy  atoms  which  may  be  the  major  mecha¬ 
nism  for  hydrogen  transport  within  the  alloy  [11,12].  As  stated  by 
Zhang  and  Angell  [13],  the  stochastic  nature  of  such  long-range 
multi-step  jumps  can  be  modeled  as  the  Fields  diffusion.  Such  a 
process  within  a  metal  particle  is  usually  simplified  as  shrinking 
core  or  contracting-envelop  model  where  the  reaction  occurs  at 
the  moving  interface  between  unreacted  core  (also  called  a-phase) 
and  reacted  shell  (also  called  (3-phase)  [  13,14].  In  the  present  macro¬ 
scopic  model,  we  deal  with  representative  volumes  which  contain 
a  sufficiently  large  number  of  metal  particles  in  each  volume.  The 
local  reaction  rate  is  the  average  over  these  particles  in  the  rep¬ 
resentative  volume.  In  addition,  the  volume  fraction  of  the  solid 
is  uniform  among  the  representative  volumes  and  hence  equal  to 
1  -  e.  The  equation  of  the  reacted  mass,  the  (3-phase,  can  then  be 
expressed  as 

(1-s)lF=Sm  (9) 

The  hydrogen  reaction  rate  is  a  function  of  temperature,  equi¬ 
librium  pressure,  and  the  portion  of  unreacted  alloy 


Sm  =  -Caexp(-|^^ln^i^(p|-/os)  (10) 

where  Ca  is  a  constant  for  specific  material  of  hydride  metals.  For 
LaNi5,  its  value  is  59.187  s-1  [9,15].  >o|  represents  the  saturated  metal 
density.  The  equilibrium  pressure  Peq  can  be  calculated  through  the 
Van’t  Floff  relationship.  Jemni  and  Nasrallah  [16]  and  Aldas  et  al.  [  17] 
adopted  the  following  form  that  directly  shows  the  dependence  of 
the  equilibrium  pressure  on  the  absolute  temperature: 

In  Peq=A-j  with./l  =  52.5  andB  =  43876.45  K  (11) 

Jemini  et  al.  [3]  further  experimentally  measured  the  P-C-T 
(pressure-concentration-temperature)  curves  at  different  temper¬ 
atures.  Starting  with  the  Van’t  Hoff  equation 


,  r»  AH 

In  Pep  =  -==- 

q  RT 


AS 

~R 


(12) 


where  the  entropy  change  AS  relates  to  the  value  of  H/M  (the 
hydrogen- to-metal  ratio),  they  determined  the  equilibrium  pres¬ 
sure  as  a  polynomial  function  of  H/M,  which  is  shown  in  Fig.  1. 
Fig.  1  indicates  that  at  low  H/Ms  Peq  decreases  rapidly  with  H/M. 


2.4.  Boundary  and  initial  conditions 


Eqs.  (1 ),  (4),  (6)  and  (9)  forms  a  complete  set  of  governing  equa¬ 
tions  with  six  unknowns:  u  (three  components),  P,  T,  and  ps.  Their 
corresponding  initial  and  boundary  conditions  are  described  as  fol¬ 
low: 

Initial  conditions: 


(13) 


Inlet  Boundaries:  The  inlet  conditions  are  set  according  to  the 
specified  temperature,  flow  rate  and  pressure. 

Walls:  No-slip  and  impermeable  velocity  condition  and  no-flux 
condition  are  applied,  while  for  the  thermal  equation  a  convection 
boundary  condition  is  adopted: 


=  /i(r-Too) 

wall 


(14) 


where  n  is  the  direction  normal  to  the  tank  outer  surface  or  the 
interface  between  the  tank  and  cooling  flow.  Note  that  from  this 


Fig.  1.  Peq  as  a  function  of  the  H/M  ratio  and  temperature.  The  numerical  fit  adopts 
the  7th  order  polynomial  suggested  by  Ref.  [19],  i.e.  a0  = -0.34863,  al  =  10.1059, 
a2  =  -14.2442,  a3  =  10.3535,  a4  = -4.20646,  a5  =  0.962371,  a6  = -0.115468, 
a7  =  0.00563776  (reproduced  from  Ref.  [19]). 


boundary  condition,  one  can  obtain  a  dimensionless  parameter  Bi 
(Biot)  number  to  compare  the  convection  at  the  tank  surface  and 
conduction  within  the  porous  bed 


(15) 


where  Lc  is  the  characteristic  length  of  the  tank.  The  heat  trans¬ 
fer  coefficient  h  is  primarily  determined  by  the  cooling  fluid.  For 
air  natural  convection  cooling,  h  is  around  l-10W(m2 1<)-1  while 
the  forced  convection  may  increase  its  value  to  30-100  W(m2  K)-1. 
For  a  small-scale  tank,  e.g.  Lc<  5  mm,  where  Bi  <  1,  temperature  is 
almost  uniform  within  the  tank.  For  larger  units  or  advanced  cooling 
systems  with  high  hs,  where  the  Bi  value  is  well  above  1,  tem¬ 
perature  varies  spatially  within  the  hydride  bed.  For  the  hydrogen 
storage  tank  considered  in  this  paper,  Bi » 1. 


3.  Numerical  procedures 


The  governing  equations  along  with  their  appropriate  boundary 
conditions  are  discretized  by  the  finite  volume  method,  with  SIM¬ 
PLE  (semi-implicit  pressure  linked  equation)  algorithm  [18].  For  the 
finite-volume  discretization,  it  is  convenient  to  unify  all  governing 
equations,  including  the  transient  terms,  in  the  following  form: 

^  +  vT(0)  =  S@  (16) 

where  0  stands  for  any  dependent  variable  in  the  governing 
equations.  Integrating  the  above  equation  throughout  an  arbitrary 
volume  V bounded  by  a  closed  surface  S,  yields: 

/  ?dt;+/  ~r{&) ' ds = j s&dv  (17) 


where  S  is  the  surface  vector.  Taking  V  and  S  to  be  the  volume  Vp 
and  discrete  faces,  Sj,  of  a  computational  cell,  respectively,  which  is 
fixed,  one  can  reach: 


/  6>di/  +  V  /  r{0)  dS=  l 

JVn  ,  JSi  Jv n 


S@dv 


(18) 
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Fig.  2.  Geometry  and  computational  domain  of  the  hydrogen  tank:  (a)  Scenario  #1 ;  (b)  Scenario  #2. 


The  final  form  of  the  discrete  finite  volume  equation  can  be 
expressed  as 

Bp0np=J2Bm@™+B(&Op)  (19) 

m 

The  above  equation  is  then  solved  by  the  algebraic  multi-grid 
(AMG)  method. 

The  geometry  and  computational  meshes  of  a  cylindrical  hydro¬ 
gen  tank  for  a  numerical  study  are  shown  in  Fig.  2.  Only  the  region 
of  z<60mm  (for  Scenario  #1)  is  filled  with  the  hydride  particles. 
The  part  of  z>60mm  is  the  expansion  region  of  the  bed.  Due  to 
the  geometrical  symmetry,  Scenario  #1  can  be  treated  as  a  two- 
dimensional  case.  To  demonstrate  the  3D  capability  of  the  model, 
we  also  consider  another  case  with  the  same  geometry  but  differ¬ 
ent  arrangement  of  the  tank,  i.e.  Scenario  #2,  or  the  tank  is  laid 
down  horizontally,  see  Fig.  2  (b).  The  masses  of  the  metal  hydride 
in  the  tanks  are  comparable  for  the  two  scenarios.  In  this  study,  we 
exclude  the  expansion  mechanism  of  alloy  particles  in  the  model 
and  assume  that  the  expansion  region  is  hollow  space  without  any 
metal  particles.  About  24,000  computational  cells  are  used  to  cap¬ 
ture  the  coupling  of  hydriding  kinetics  and  heat/mass  transport. 
Geometrical  and  operating  parameters  of  the  storage  tank  as  well 
as  physical  properties  are  listed  in  Table  1.  Adaptive  time  stepping 


Table  1 

Geometrical,  physical,  and  operating  parameters. 


Quantity 

Value 

Tank  dimension,  height/diameter 

80/50  mm 

Initial/inlet  temperature,  T0/Tin 

20/20  °C 

Inlet  pressures,  P 

10.0  bar 

Porosity/permeability  of  the  hydride  bed,  e 

0.3/10-8  m2 

Tortuosity  of  the  solid  matrix  and  gas,  rg/rs 

1.5/1. 5 

Thermal  conductivity  of  the  bed,  kglks 

0.167/2.0  Wm-1  K-1 

Specific  heat  capacities,  c^/csp 

14.89/0.419  kj  kg-1  K_1 

Activation  energy,  Ea 

21179.6  J  mol”1 

Kinematic  viscosity  of  hydrogen,  vg 

1.03  x  10“4  m2  s-1 

Density  of  the  hydride,  ps 

4200  kg  m-3 

Cooling  fluid  temperature,  T » 

20°C 

Heat  transfer  coefficient,  h 

600  Wm-2  K_1 

is  used  in  which  the  current  time  step  is  inversely  proportional  to 
the  magnitude  of  the  maximum  reaction  rate  at  the  previous  time 
step  with  the  maximum  of  0.1  s.  For  iterations  in  each  time  step, 
equation  residuals  are  all  smaller  than  10-7.  A  typical  simulation 
takes  ~1  h. 

4.  Results  and  discussion 

Fig.  3  compares  the  model  prediction  and  experimental  data. 
Detail  of  the  experimental  setup  can  be  found  in  Ref.  [3  ].  The  reactor 
tank  is  cylindrical  with  geometrical  parameters  and  alloy  particle 
arrangement  the  same  as  Scenario  #1.  The  tank  was  submerged  in 
a  bath  at  a  constant  temperature  (Too).  In  their  experiment,  three 


t(s) 


Fig.  3.  Validation  with  experimental  data  in  terms  of  temperature  evolutions  (Sce¬ 
nario  #1).  Experimental  data  are  from  Ref.  [3].  Point  A:  located  at  r=15mm  and 
z  =  35  mm;  Point  B:  located  at  r=  15  mm  and  z  =  25  mm;  Point  C:  located  at  r=  15  mm 
and  z=  15  mm. 
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Fig.  4.  The  temperature  distributions  in  a  cross-section  of  the  tank  at  1  s,  5  s  and  10  s  for  Scenario  #1. 


thermal  couples  were  placed  within  the  alloy  bed  (location  A,  B,  and 
C)  to  measure  the  local  temperature  evolution.  This  figure  shows  a 
satisfactory  agreement.  At  the  beginning,  the  model  over-predicts 
the  temperatures,  which  may  be  due  to  neglecting  the  heat  capac¬ 
ity  of  the  tank  wall  in  the  simulation.  The  validation  is  based  on 
Scenario  #1,  where  hydrogen  is  injected  from  the  above. 

Figs.  4  and  5  present  the  evolution  of  temperature  distributions 
within  the  tank  for  Scenario  #1.  Fig.  4  shows  the  temperature  dis¬ 
tributions  at  the  beginning  of  the  fueling  (1,  5  and  10  s).  It  can  be 
seen  that  the  temperature  in  the  porous  hydride  bed  (z<  60  mm)  is 
almost  uniform  and  increases  rapidly  from  the  initial  temperature. 
At  1  s,  the  upper  corners  near  the  interface  between  the  alloy  bed 
and  expansion  space  shows  different  gradient  patterns.  This  is  due 
to  the  impact  of  hydrogen  flow  cooling,  which  will  be  explained  in 
detail  by  the  flow  fields  plotted  later.  At  5  s,  it  can  be  observed  that 
there  exist  obvious  temperature  gradients  near  the  wall  and  upper 
surface  of  the  metal  bed.  In  particular  near  the  lower  corners,  the 
cooling  effect  cannot  be  neglected.  As  the  hydride  process  proceeds, 
the  non-uniformity  of  the  temperature  distribution  becomes  evi¬ 
dent  in  the  bed.  At  100  s,  the  area  near  the  wall  is  much  colder  than 
the  core  region.  This  is  partly  due  to  the  slow-down  overall  reaction 
rate  arising  from  the  depletion  of  the  reactant  metal,  which  reduces 


the  heat  production  rate  and  hence  wall  temperature.  At  2000  s,  the 
high  temperature  only  concentrates  in  the  small  region  of  the  core. 

Figs.  6  and  7  present  the  FI/M  ratio  distributions  at  the  six  time 
instants.  Fig.  6  shows  that  the  H/M  ratio  is  almost  uniform  at  1  s, 
while  higher  values  appear  near  the  wall  and  the  bed  upper  surface 
at  5  and  10  s  but  the  spatial  variation  is  fairly  small  in  the  hydride 
bed.  Fig.  7  shows  the  high  H/M  region  increases  its  area  towards  the 
core,  where  remains  a  low  H/M  ratio  in  most  of  the  time.  Even  at 
2000  s,  the  hydrogen  content  in  the  alloy  near  the  core  just  reaches 
about  half  of  the  saturation  value. 

Figs.  8  and  9  display  the  flow  fields  at  six  time  instants.  It  can 
be  seen  that  the  initial  velocities  are  much  higher  than  the  ones 
at  the  latter  stages.  At  the  first  second,  the  velocity  is  able  to  reach 
0.2  ms-1,  indicative  of  a  fast  reaction  rate  at  the  initial  stage.  The 
fast  kinetics  is  benefited  by  the  initial  lower  temperature  and  H/M 
ratio.  Note  that  the  low  H/M  ratio  also  affects  the  equilibrium  pres¬ 
sure,  Peq  (see  Fig.  1 ).  In  the  expansion  space,  the  velocity  magnitude 
remains  nearly  constant,  while  the  velocity  in  the  bed  decreases 
when  moving  downwards  due  to  the  absorption.  Similar  trends  are 
also  shown  at  other  time  instants.  At  the  very  end  stage,  i.e.  2000  s, 
there  seems  a  hydrogen  mass  sink  in  the  core  region,  which  can 
be  explained  by  the  local  higher  content  of  remaining  unsaturated 
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Fig.  5.  The  temperature  distributions  in  a  cross-section  of  the  tank  at  100  s,  1000  s  and  2000  s  for  Scenario  #1. 
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Fig.  6.  The  H/M  ratio  or  (H/M)/(H/M)S  distribution  in  a  cross-section  of  the  tank  at  1  s,  5  s  and  10  s  for  Scenario  #1. 
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Fig.  7.  The  H/M  ratio  or  (H/M)/(H/M)S  distributions  in  a  cross-section  of  the  tank  at  100  s,  1000  s  and  2000  s  for  Scenario  #1. 
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Fig.  8.  The  velocity  distributions  in  a  cross-section  of  the  tank  at  1  s,  5  s  and  10  s  for  Scenario  #1. 
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Fig.  9.  The  velocity  distributions  in  a  cross-section  of  the  tank  at  100  s,  1000  s  and  2000  s  for  Scenario  #1. 


metal,  as  shown  in  Fig.  7.  In  addition,  except  the  lsts,  the  flow 
profiles  in  the  expansion  region  clearly  indicate  the  boundary  layer 
(BL)  development  near  the  wall.  This  is  due  to  the  fact  that  the 
BL  thickness  is  linear  with  the  reverse  of  square  root  of  Re  (the 
Reynolds  number).  At  the  1  st  s,  Re  is  relatively  high  due  to  the  large 
flow  rate,  leading  to  an  almost  invisible  BL,  while  in  other  cases  par¬ 
ticular  at  latter  stage  when  the  velocity  is  low,  the  velocity  profile 
near  the  wall  in  the  expansion  region  is  evident.  In  the  alloy  bed, 
the  Darcy’s  force  dominates  and  the  wall  effect  becomes  negligible. 
Therefore,  when  entering  the  bed,  part  of  the  hydrogen  mass  out  of 
the  boundary  layer  goes  towards  the  wall  for  the  local  adsorption,  as 
indicated  in  these  two  figures.  This  lateral  flow  provides  extra  cool¬ 
ing.  For  the  strong  flow  case  at  the  1  st  s,  the  cooling  effect  is  evident, 
recasting  the  local  temperature  distribution  at  the  upper  corners 
of  the  hydride  bed  as  indicated  by  the  first  plot  in  Fig.  4.  Further,  a 
dimensionless  number,  the  Peclet  number  (=Re  x  Pr  =  uR/a  where  a 
is  the  thermal  diffusivity),  can  be  used  to  evaluate  the  importance  of 
the  cooling  effect  provided  by  the  hydrogen  flow.  For  u  ~  0.1  m  s-1 , 
R~  25  mm,  and  1.0  x  10-4  m2  s-1,  the  Peclet  number  is  ~25, 
thus  convection  is  dominant.  Even  with  the  velocity  decreasing 
to  ~0.01  ms-1  (e.g.  at  100s),  convection  effect  is  still  comparable 


with  the  conduction  (the  Peclet  number  ~2.5).  At  1000  s,  the  Peclet 
number  is  ~0.5  in  the  bed  near  the  inlet,  indicating  the  convection 
plays  an  important  role  in  the  local  heat  transfer. 

Fig.  10  plots  the  temperature  evolutions  of  several  typical  points 
in  the  hydride  bed.  It  can  be  seen  that  the  temperatures  near  the 
wall  drop  fast  in  the  initial  stage  (Note  that  there  is  a  jump  in  tem¬ 
perature  at  the  very  beginning  when  the  exothermal  reaction  heats 
up  the  bed  as  shown  in  Fig.  4.  Flowever,  this  part  is  relatively  short 
and  not  clearly  visible  in  this  figure.  We  exclude  our  discussion  on 
this  period  in  this  paragraph).  In  particular,  the  one  at  the  bottom 
corner  decreases  fastest.  After  the  initial  fast  drop,  the  decreas¬ 
ing  pace  slows  down.  In  contrast,  the  two  locations  in  the  core 
region  indicate  a  totally  different  trend.  Initially,  their  temperatures 
remain  high  with  little  change,  followed  by  sharp  drops  after  ~1500 
and  1900  s,  respectively.  Fig.  11  shows  the  evolutions  of  the  mass 
absorption  rates  at  these  locations.  It  can  be  seen  that  the  ones  near 
the  wall  are  much  faster  than  the  ones  near  the  core  in  the  initial 
100  s  and  decrease  nearly  exponentially  with  time.  As  to  the  two 
locations  near  the  core,  similar  trends  are  indicated:  first  decrease 
partly  due  to  the  depletion  of  the  metal,  then  increase  upon  the 
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Fig.  12.  The  temperature  distributions  in  the  tank  at  four  time  instants  for  Scenario  #2. 


temperature  drop,  followed  by  another  reduction  caused  by  the  fact 
that  the  metal  depletion  dominates  again.  In  addition,  the  middle 
one,  i.e.  z  =  30  mm  and  r = 0  mm,  falls  sharply  at  the  very  beginning, 
several-order-of-magnitude  lower  than  others,  which  is  due  to  the 
local  high  temperature.  In  addition,  even  though  the  reaction  at  the 
other  point,  i.e.  z  =  60  mm  and  r = 0  mm,  is  fast  in  most  time,  the  tem¬ 
perature  is  lower  as  shown  in  previous  figure,  which  is  benefited 
by  the  cooling  of  the  hydrogen  flow. 

As  mentioned  before,  Scenario  #1  is  essentially  a  two- 
dimensional  problem,  while  the  problem  of  Scenario  #2  is 
three-dimensional.  Fig.  12  shows  the  3D  temperature  distributions 
for  Scenario  #2  at  four  time  instants.  Similar  to  the  former  case, 
the  bed  temperature  is  almost  uniform  initially,  while  varies  spa¬ 
tially  later  with  the  core  at  the  highest  temperature.  Fig.  13  shows 
the  velocity  distributions.  Different  with  Scenario  #1,  the  velocity  in 


the  expansion  is  much  higher  than  the  one  in  the  hydride  bed.  In  the 
expansion  space,  the  velocity  decreases  along  the  flow  due  to  the 
bed  absorption.  In  the  bed,  the  flow  direction  indicates  that  hydro¬ 
gen  is  supplied  from  both  the  inlet  and  expansion  region.  At  2000  s, 
hydrogen  flow  is  drawn  towards  the  core  where  the  reaction  is  still 
active,  similar  to  Fig.  9.  In  addition,  high  velocity  in  the  hydride 
bed  provides  extra  cooling  for  the  exothermic  reaction  as  shown 
in  Scenario  #1.  However,  repeated  hydrogenation/dehydrogenation 
cycling  will  cause  metal  hydride  alloy  particles  pulverized.  The 
hydrogen  flow  or  other  external  forces  such  as  mechanic  shock  or 
vibration  may  bring  the  resultant  fine  particles  to  the  bottom  of  the 
cylinder  forming  a  denser  layer.  This  layer  will  have  much  lower  per¬ 
meability  due  to  the  small  pore  dimension  (see  Eq.  (3)),  affecting 
the  hydriding  process.  These  two  effects  must  be  balanced  through 
sophisticate  design  in  practice. 
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Fig.  13.  The  velocity  distributions  in  the  tank  at  four  time  instants  for  Scenario  #2:  (a)  the  middle  cross-section  of  the  tank;  (b)  the  middle  cross-section  of  the  hydride  bed. 


5.  Conclusions 

A  three-dimensional  model  was  developed  to  investigate  the 
coupled  mass/heat  transport  phenomena  during  the  absorption 
process  in  hydrogen  storage  tanks.  The  model  consists  of  principles 
of  mass,  momentum,  and  energy  conservations  and  the  absorption 
kinetics  of  the  metal  hydride.  Numerical  simulations  were  carried 
out  to  reveal  the  spatial  distributions  and  evolutions  of  the  key 
quantities  within  a  cylindrical  tank  for  two  scenarios.  The  results 
indicated  that  hydrogen  absorption  is  almost  uniform  at  the  ini¬ 
tial  stage  but  varies  spatially  at  the  latter  one.  In  the  core  region, 
temperature  remains  high  for  quite  a  long  time  due  to  insufficient 
cooling,  followed  by  a  fast  drop  before  the  absorption  is  completed. 
The  initial  absorption  is  fast  due  to  the  initial  low  temperature  and 
H/M  ratio.  This  fast  reaction  results  in  a  hydrogen  flow  in  the  bed 
with  the  velocity  magnitude  much  higher  than  the  latter  stage.  The 
strong  flow  also  provides  extra  cooling  to  the  bed  near  the  inlet. 
Further,  we  found  that  the  temperature  increases  rapidly  and  uni¬ 
formly  at  the  beginning.  It  also  implies  that  an  optimal  thermal 


design  is  important  for  achieving  fast  refueling  of  hydride  storage 
tanks.  Future  work  includes  investigating  the  impacts  of  parame¬ 
ters  on  hydrogen  fueling,  modeling  hydrogen  desorption  process, 
and  integrating  the  hydrogen  tank  model  with  our  fuel  cell  ones 
[20,21]. 
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